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ABSTRACT 

We create realistic, full-sky, half-arcminute resolution simulations of the microwave sky matched to the most 
recent astrophysical observations. The primary purpose of these simulations is to test the data reduction pipeline 
for the Atacama Cosmology Telescope (ACT) experiment; however, we have widened the frequency coverage 
beyond the ACT bands and utilized the easily accessible HEALPix map format to make these simulations 
applicable to other current and near future microwave background experiments. Some of the novel features 
of these simulations are that the radio and infrared galaxy populations are correlated with the galaxy cluster 
and group populations, the primordial microwave background is lensed by the dark matter structure in the 
simulation via a ray-tracing code, the contribution to the thermal and kinetic Sunyaev-Zel'dovich (SZ) signals 
from galaxy clusters, groups, and the intergalactic medium has been included, and the gas prescription to model 
the SZ signals has been refined to match the most recent X-ray observations. The cosmology adopted in these 
simulations is also consistent with the WMAP 5-year parameter measurements. From these simulations we 
find a slope for the 1200-^^200 relation that is only slightly steeper than self-similar, with an intrinsic scatter 
in the relation of ^ 14%. Regarding the contamination of cluster SZ flux by radio galaxies, we find for 148 
GHz (90 GHz) only 3% (4%) of halos have their SZ decrements contaminated at a level of 20% or more. We 
find the contamination levels higher for infrared galaxies. However, at 90 GHz, less than 20% of clusters with 
M200 > 2.5 X lO'^^M© and z < 1.2 have their SZ decrements filled in at a level of 20% or more. At 148 GHz, 
less than 20% of clusters with M200 > 2.5 x IO'^^Mq and z < 0.8 have their SZ decrements filled in at a level of 
50% or larger. Our models also suggest that a population of very high flux infrared galaxies, which are likely 
lensed sources, contribute most to the SZ contamination of very massive clusters at 90 and 148 GHz. These 
simulations are publicly available and should serve as a useful tool for microwave surveys to cross-check SZ 
cluster detection, power spectrum, and cross-correlation analyses. 

Subject headings: cosmic microwave background - galaxies: clusters: general - galaxies: general - intergalac- 
tic medium - large-scale structure of universe - methods: simulations 



1. INTRODUCTION 

Simulations of the microwave sky play a crucial role in 
accurately extracting constraints on cosmology potentially 
measurable by current and upcoming microwave background 
observatories, such as the Atacama Cosmology Telescope 
(ACT)^ the South Pole Telescope (SPT)'', and the Planck 
satellite'". These surveys will have complex data process- 
ing pipelines and analysis strategies, and any artifacts caused 
by any one of the many analysis steps can best be controlled 
by simultaneously running simulations of the microwave sky, 
with known cosmological and astrophysical inputs, through 
the same analysis pipeline. These microwave background ex- 
periments, which will have finer resolution and higher sensi- 
tivity than achieved to date, can potentially constrain cosmo- 
logical parameters via their measured microwave power spec- 
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trum, higher order correlation functions, Sunyaev-Zel'dovich 
(SZ) cluster detections, and cross-correlation studies of vari- 
ous microwave components. Moreover, these microwave data 
sets can be combined with surveys at other wavelengths to 
yield additional constraints. As such, the accurate recovery 
of known cosmological inputs, through the variety of meth- 
ods mentioned above, must be tested with simulations. This 
necessitates high-resolution, cosmological-scale, microwave 
simulations that include the most accurate astrophysical com- 
ponents as our current observations will allow. 

The microwave sky consists primarily of the following 
components: (1) the primordial microwave background, 
which is lensed by the intervening structure between the last 
scattering surface and observers today, (2) the thermal SZ sig- 
nal from the hot gas in galaxy clusters, groups, and the inter- 
galactic medium, (3) the kinetic SZ signal from the bulk mo- 
tion of galaxy clusters, groups, and the intergalactic medium, 
(4) higher-order relativistic corrections to the thermal and ki- 
netic SZ signals, (5) a population of dusty star forming galax- 
ies that emit strongly at infrared wavelengths but have signif- 
icant microwave emission, (6) a population of galaxies that 
emit strongly at radio wavelengths, including active galactic 
nuclei, which also have significant emission at microwave fre- 
quencies, and (7) dust, synchrotron, and free-free emission 
from the Galaxy. 

The simulations discussed here include the above compo- 
nents, mapped at six different frequencies: 30, 90, 148, 219, 
277, and 350 GHz. The frequencies 148, 219, and 277 ai-e the 
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central ACT observing frequencies, as these simulations are 
used by the ACT team to inform their data analysis pipeline. 
The additional frequencies have been added to widen the 
applicability of these simulations to other microwave back- 
ground experiments. The only components listed above that 
have been excluded from the microwave maps are the Galac- 
tic synchrotron and free-free emission. These two compo- 
nents are not expected to be significant at ACT frequencies, 
and we leave their inclusion to future work. The maps them- 
selves are full-sky HEALPix maps of resolution 0.4 arcmin- 
utes (Nside=8192) and consist of full-sky maps of each com- 
ponent separately as well as a combined sky map, at each 
frequency. These HEALPix maps can be interpolated into 
flat sky maps using the bilinear interp olation subroutin e in- 
cluded in the HEALPix distribution CGorski et al.ll2005h . In 
addition, catalogs of all the halos in this simulation and their 
various properties, as well as catalogs of the individual in- 
frared and radio galaxies that populate the halos, are provided. 
These simulations are public and can be downloaded from 
|http://lambda. gsfc.nasa.gov/toolbox/tb_cmbsim_ov.cfm 

There are of course ma ny simulations of components of 
the microwave sky (e.g. Gawiser & Smoot 1997; Metzler' 
1998; Gawiser etal. 1998; Sokasian et al. 2001; White et aL 
~002; daSilvaetal. 2004; Diaferio et al. 2005; Mo tl et ai] 



2005; 
2007; 



Nagai 2 006; .Hallman et a l. 2007 



Roncare lli et al. 120071: iMalte Schafer & Bartelmann 



2007c iLeach et alj|200 8!: 'Pace et al. 2008; Shaw et al 2QQ8t 



Wik et al.ii2008t iCarb one et al. 2009; Peel et al. 2009|. The 



simulations discussed here include the components men- 
tioned above, in an internally consistent manner, matched 
to the most recent observations. These simulations have 
improved upon t he previous simulations discussed 
ISeh^al et alj (l200'7h in the following ways: 



m 



• radio and infrared galaxy populations are clustered and 
correlated with the galaxy cluster populations 

• the primordial microwave background has been lensed 
by the dark matter structure in the N-body simulation 
via a ray-tracing code 

• the contribution to the kinetic SZ signal from the inter- 
galactic medium has been included, in addition to that 
from galaxy clusters 

• the gas prescription to model the SZ signals has been 
refined to match the most recent X-ray observations 

• maps have been created at a wider range of frequencies 
to increase the applicability to other microwave back- 
ground experiments 

• standard HEALPix format for the maps has been 
adopted 

We give a brief summary of the key elements of these 
simulations below, which we expand upon in greater detail 
throughout §2 and §3. 

The large-scale structur e simulation was carried out using 
a tree-particle-mesh code dBode et alj|2000l: iBode & Ostrikerl 
l2003h . with a simulation volume of 1000 h~^Mpc on a 
side containing 1024^ particles. The cosmology adopted 

is {nh,n,n,nA,h,n,, (Tfi) = (0.044, 0.264, 0.736, 0.71, 

0.96, 0.80), consistent with the WMAP 5-year results 



(iKomatsu et al.ll2009h . The mass distribution covering one 
octant of the full sky was saved, and halos with a friends- 
of-friends mass above 1 x IO^^Mq and with a redshift below 
z = 3 are identified. The octant is then replicated to create 
a full-sky simulation. Thus the large-scale structure compo- 
nents in these simulations are unique for one octant, while the 
primary microwave background and Galactic dust are mod- 
eled for the full sky. 

The primary microwave background has the same cosmol- 
ogy as the large-scale structure and matches the WMAP 5- 
year Internal Linear Combination map (Bennett et al. 20031 
Hinshaw et al. 2009) on large-scales (l < 20). The microwave 
background is then lensed by the l arge-scale structure in the 
simulation with a ray-tracing code dPas & Bodell2008h . 

The Sunyaev-Zel'dovich signal is derived by adding a gas 
prescription to the N-body halos which assumes a polytropic 
equation of state and hydrostatic eq uilibrium. T his model, 
which is described in more detail in iBode et aP ( 12009 ). ad- 
justs four free parameters which are calibrated against X-ray 
gas fractio ns as a function of temperature from the sample of 
ISun et all (j2009) and Vikhlinin et al. (2006). These four free 
parameters are the star formation rate (which depends on the 
halo mass and redshift), non-thermal pressure support (from 
cosmic rays or magnetic fields), dynamical energy transfer 
from dark matter to gas (a depletion factor which lowers the 
average baryon fraction to less than the cosmic mean), and 
feedback from active galactic nuclei (which is tied to the star 
formation rate). This model describes the gas prescription for 
halos with friends-of-friends mass > 3 x lO'^M© and z < 3. 

Smaller halos down to 1 x 1O'^M0 and the intergalactic 
medium for z < 3 are treated with a different method. The 
velocity dispersion of a given simulation particle is calculated 
using its neighbors, and this velocity dispersion is related to 
the pressure (and thus the thermal SZ signal) via a propor- 
tionality co nstant derived from hydrodynamical simulations 
(iTrac et a l.l l2b09i) . The kinetic SZ is calculated from the line- 
of-sight momentum of the particles. For the 3 < z < 10 Uni- 
verse, the temperature is assumed to be isothermal at 10* K, 
and mass density and momentum shells at various redshift 
slices are used to construct the thermal and kinetic SZ signals. 

The infrared source model has six parameters used t o pop- 
ulate halos, and is partially based on Righi et al. ("2008). This 
model satisfies the following observational constraints: the 
COBE/FIRAS background as in lFixsen et al.i ( 1998) is an up- 
per limit on the total infrared intensity, the effective spec- 
tral index for the spectral intensity is ^ 2.6 be tween 145 
and 350 GHz dKnox et all 120041: iDunne et al.ll20bQ) . source 
counts are compat ible with SCUBA counts at 353 GHz 
('Copp in et alJl2006l) . source counts and clustering are com- 
patible with BLA ST me asured cou nts and power spectrum 
(iPatanchon et all "2009: V iero et all [2009). and the power of 
the infrared sources is below the point source upper limit de- 
rived by ACBA R at 150 GHz fixing tig to the WMAP 5-yeai- 
measurement (iReichardt et alj|2009l) . 

The model for the radio sources is made by creating a halo 
occupation distribution and luminosity distribution for radio 
sources at 151 MHz to match the 151 MHz radio luminos- 
ity function (RLE) at z ~ 0.1. 151 MHz was chosen because 
at such low frequencies the RLE is dominated by the steep 
spectrum lobes of the radio sources (as opposed to the flat- 
ter spectrum cores) and is thus less sensitive to biases due to 
orientation effects. Once the 151 MHz RLE is matched at 
low redshifts, the radio sources are divided into two popu- 
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lations, similar to the types I and II of the iFanaroff & RilevI 
( 119741) classification. These are given differing number den- 
sity evolutions with redshift, and the RLF at higher redshift is 
matched to observations. We use a relativistic beaming model 
to separate how much of a source's flux at 151 MHz is from 
the core as opposed to the lobes, and we adjust the shape of the 
radio core spectral energy distribution to mat ch source count s 
at J/ » 151 MHz, where the core dominates tLin et alj2009al) . 

The Galactic dust emission maps are from t he "model 8" 
prediction from lFinkbeiner et alJ (Il999h . which'Benne tt et aT] 
P003) show to be a reasonable template for dust emission in 
the WMAP maps. 

In ^ we discuss the construction of each component of 
these simulations and comparison to observations in more de- 
tail. In ^ we discuss further properties of the simulated 
maps, and in ^ we summarize and conclude. 

2. SIMULATIONS 

2. 1 . Primary Unlensed CMB 

We follow the method of ISehgal et alJ ([2003), construct- 
ing a primary microwave background map with the observed 
large-scale structure and small scale structure consistent with 
theoretical expectations. For large scales, we u se the WMAP 
5-year Internal Linear Co mbination (ILC) map jBennett et alj 
12003; Hins haw et al.1 '2009). For / < 20, we take the harmonic 
coefficients (a;,,,) from the ILC map with no modification. At 
smaller scales, where the ILC map is smoothed significantly, 
a Gaussian random realization is added, so that the ensemble 
average power equals the theoretical power spectrum which 
maximizes the WMAP 5-year likelihood. Harmonic coeffi- 
cients are computed to /^ax = 8192, and we synthesize a map 
with HEALPix parameter A^^ide = 4096 (0.859' pixels). Then 
we convert to a map with HEALPix parameter N^ias = 8192 
(0.429' pixels) using the a 1 1 e r_a Im subrou tine included in 
the HEALPix distribution ( iGorski et alj2005l) . The lower res- 
olution (A^side = 4096) map is us ed to generate the lensed mi- 
crowave background map (see § l2.3l l. as that is more efficient 
than working with the A^side = 8 192 map directly and results in 
no significant loss of accuracy. The lensed map is then also 
converted to an A^side = 8192 map. 

2.2. Large-Scale Structure Simulation 

A simulation of the large-scale structure of the Universe is 
performed usi ng the same methods as previo usly described 
in iBode et all (120071) and ISehgal et all (120071) . and we refer 
the reader to their papers for additional details. Within a pe- 
riodic box of comoving side length L = 1000 /i"'Mpc, 1024^ 
dark matter particles a re evolved using the Tree -Particle-Mesh 
(TPM) N-body code dBode & Ostrikeri 12003b . The paiticle 
mass is 6.82 x 10'" h~^MQ and the gravitational spline soft- 
ening length is set to e = 16.276 /?~'kpc. Th is simulation has 
been previously used by Bod e et aH (l2009h in developing a 
model for the intracluster medium. 

A putative observer at redshift zero is placed at one corner 
of the periodic box, and the mass distribution along a past 
light cone, covering one octant (box coordinates x,y,z > 0) 
of the sky, is saved as the simulation is running. To create 
full-sky maps this octant is replicated, as described in |3] 

At each simulation time step, spanning the redshift range 
Zi > z > Z2, particles are found in a thin spherical shell at the 
comoving distance d corresponding to the light travel time to 
this observer, spanning the range d(zi) > d > dizi). For z < 3 
all of these particles are saved to disk and for z < 10 pixelized 
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Fig. 1 . — The halo abundance above a minimum mass per unit redshift per 
square degree measured from the large-scale structure simulation (p oints) is 
compared with the semi-analytic fitting formula (lines) from Jenki ns et alj 
120011) . From top to bottom, the four sets of points are for the minimum 
masses M„i„ = 6.82 x lO'-, 2 x lO'^ 1 x lO", and 3 x lO'* /r'MQ. 

versions of the shells are saved, as described below. There 
are 579 such shells in all; the range Az = Zi -Z2 of each shell 
depends on the time step size set by the TPM code. It de- 
creases from Az ~ 0.08 at z ~ 9 to Az ~ 0.03 at z ~ 3, and to 
Az < 0.005 for z < 0.8. For comoving distances larger than 
the box size, d > 1000 /z"'Mpc, there can be some duplica- 
tion in structures as the periodic box is repeatedly tiled. This 
tiling is done without random rotations to preserve the period- 
icity and prevent discontinuities in the large-scale structure. 
However, any repetition is usually at a different redshift and 
evolutionary state. In cases where a halo or filament appears 
twice at similar redshifts, it is viewed at two different angles, 
making each projected image unique. 

For all the shells (z < 10) particles are subdivided by angu- 
lar coordinates using the HEALPix scheme for pixelation of 
a sphere. The projected mass and momentum in each pixel of 
the shell are com puted and saved to disk (as described previ- 
ously in lDas & B ode 2008). These pix elized shells are used 
to model the gravitational lensing 02.3I ) and high-redshift SZ 
distortion 02.4l i of the microwave background. 

For lower redshifts z < 3, the positions and velocities of 
all particles in the light cone are also saved to disk for post- 
processing. The approximately 55 billion particles thus saved 
occupy 1.3 TB of disk space and are used to model the SZ 
effect ( ^2.4b in detail. A friends-of-friends (FoF) halo finder, 
with a standard linking length /jfoF equal to 1/5 of the mean 
interparticle spacing, is used to identify dark matter halos and 
tag the particles belonging to them. This halo catalog is used 
to model infrared 02.5t and radio 02.6I ) point sources. 

As shown in Fig.[Tl the halo mass functions measured from 
particles in the ligh t cone agree well wi th the semi-analytic 
fitting formula from I Jenkins et al.l (1200 Ih for a linking length 
bpoF = 0.2. The number of halos above a minimum mass 
per unit redshift per square degree is plotted for four min- 
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imum masses: M^^^ = 6.82 x 10'^ 2 x 10'^ 1 x 10'^ and 
3 X 10'"^ hr^MQ. The first limit, the equivalent of 100 parti- 
cles, is the lowest mass of halos that are populated with point 
sources. The second value corresponds to massive halos with 
enough resolution to model the intrahalo gas in greater detail. 
For any cluster survey, the minimum detectable mass is likely 
to be a function of redshift, and the last two minimum mass 
values reflect the range expected for upcoming SZ surveys. 
The data points in Fig. [T] are measured using all halos in the 
octant to minimize sample variance, but then are normaUzed 
per square degree. 

2.3. Lensing of the Primary CMB 

Gravitational lensing is a secondary anisotropy of the mi- 
crowave background caused by the deflection of primary 
microwave background photons by intervening la rge-scale- 
structure potentials (see iLewis & Challinoii l2006l for a re- 
cent review). Lensing produces non-Gaussianities, smooths 
out features in the power spectrum and moves power from 
large to small scales. High resolution microwave background 
experiments, such as ACT, SPT, and Planck, will probe the 
scales corresponding to the multipole range 500 < ^ < 10000 
with unprecedented precision. It will be necessary to in- 
clude the effect of lensing for performing precision cosmo- 
logical tests with such datasets. By studying the lensing dis- 
tortions in the microwave sky, we should also be able to re- 
construct the projected matter distribution or the de flection 
field in the Universe (Hu & Okamoto 2002; Okamot o & Hul 
|2003; Hirata & Seljak 2003; Yoo & Zaldarriaga 2008). Such 
reconstruction will help break parameter degeneracies in the 
microwave background, maki ng it sensitive to dark energy dy- 
namics and neutrino mass (iSmith et al.ll2008t Ide Putter et alj 
|2009). Also, the cross-correlation of microwave background 
lensing with various tracers of large scale structure will lead to 
better understanding of galaxy environments, structure forma- 
tion, geometry and gravity on large scales (lAcauaviva et alJ 
I2OOI iDas & SpergelllMiQt iValUnotto et alJl2009h . A realis- 
tic simulation of microwave lensing internally consistent with 
(i.e. based on the same large-scale-structure simulation as) the 
other simulated secondary components, such as the thermal 
and kinetic SZ effects and infrared and radio point sources, is 
essential for training the power spectrum and parameter esti- 
mation pipeline. Such simulations are also necessary ingre- 
dients in developing realistic lens reconstruction algorithms 
dAmblard et al. 2004; Perotto et al. 2009). Here we present 
the salient aspects of a full-sky microwave background lens- 
ing simulation at arcminute resolution. 

A simulation of the lensed microwave background es- 
sentially consists of remapping the points of the primary 
microwave background sky according to deflection angles 
derived from a simulated convergence field. Several ap- 
proaches have been proposed for implementing the remap- 
ping as well as for generating the de flection field (Lewis 20051 
Das & Bodd 120081 ICarbone et aP 120071 I2OO8; Basak et al] 



2008h . In this work, we foUow Das & Bode (200|) in aU re- 



spects, except that we create a full-sky convergence map (in- 
stead of the polar cap geometry adopted in that paper) and 
lens a full-sky primary microwave bac kground map . The d e- 
tails of the algorithm can be found in iDas & Bodd (120081) - 
here we delineate the main steps: 



Simulation 
Theory 
Linear Theory 




1. As described in ^2.21 at each time step of the large- 
scale-structure simulation, the positions of the particles 
in a thin shell (of width corresponding to the time step) 



Fig. 2. — Power spectnim of the simulated effective kappa map (dots) com- 
pared to tlieory. with (solid line) and without (dashed line) nonlinear correc- 
tions to the matter power spectrum. 

are recorded onto the octant of a HEALPix pixelated 
sphere. Knowing the particle mass and the number of 
particles falling into each pixel, we generate a surface 
mass density plane (expressed in mass per steradian) by 
dividing the total mass in a pixel by the solid angle sub- 
tended by it. We then subtract the cosmological mean 
of the surface density to generate a surface overdensity 
on each shell. 

2. We then generate the effective convergence map (k) 
up to the maximum redshift in the large-scale-structure 
simulation (z= 10) by combining the surface overdensi- 
ties in the octant shells with pr oper geometrical w eight- 
ing (see, eqns. 15 and 20 in iDas & Bod^ 12008) . To 
accurately lens the microwave background, we need to 
include the contribution to the convergence from red- 
shifts beyond the farthest simulated shell (z > 10). We 
do this by adding in a Gaussian random realization for 
this extra convergence based on a theoretical conver- 
gence power spectrum (calculating the contribution to 
the convergence power from z > 10 structure). At the 
end of this process, we have a realization of the effec- 
tive convergence K{n) on the octant. 

3. We replicate the octant to create a full-sky convergence 
map. The detailed mapping of the octant to the full-sky 
is described in ^ 

4. We invert the relation, k= jV^cj), between the conver- 
gence and the effective lensing potential (p in harmonic 
space to obtain the harmonic components of 0, and 
generate the deflection field, a.(n), using the relation 
a =-V0. 

5. We use the deflection field to find the source-plane po- 
sitions, 6s, corresponding to the observed positions, 9, 
of the rays (centers of pixels) according to: 6s = 6+0.. 

6. Finally, we sample the unlensed microwave back- 
ground at the source-plane positions using an accu- 
rate interpolation technique to produce the lensed mi- 
crowave background map. 
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Fig. 3. — Left: Lensed and unlensed microwave background power spectra obtained from the simulation compared with the theoretical models. The blue dots 
represent the power spectrum obtained from the simulated lensed microwave background map while the red solid line represents the theoretical lensed microwave 
background power spectrum from CAMB. The red dots represent the power spectrum of the input unlensed map and the green soHd Hne is the theoretical input 
power spectrum. Note that the simulation and theory lines coincide to high precision. Right: A smaller section of the spectra on the left. Note that the spectrum 
is multiplied by the fourth power of the multipole to make clearer the difference between lensed and unlensed spectra. 




R.A. 02000) 

Fig. 4. — Difference between the lensed and the unlensed primary mi- 
crowave background map in a ~ 10 degree X 10 degree patch. 



For the results presented here, we have performed the above 
steps at the HEALPix resolution A^side = 4096. Then we con- 
vert to a map with HEALPix parameter N^^^ = 8 192 using the 
alt er_alm subroutin e included in the HEALPix distribu- 
tion (lG6rski et al.ll2005h . 

Fig. 12] compares the power spectrum of the simulated ef- 
fective K map with two theoretical models: o ne with and 
the ot her without nonlinear corrections (following Smith et al.l 
I2OO3I) to the matter power spectrum. The theoretical spec- 
tra have been computed with the publicly available CAMB" 
code. As expected, the simulated power spectrum agrees bet- 
ter with the theory curve with nonlinear corrections. There 
are two things worth noticing in the plot: the large scatter at 

" www.camb.info 



the low multipole end of the power spectrum, and the excess 
power in the simulation beyond t ^ 2000. The scatter is ex- 
pected because we generated the full-sky convergence map by 
replicating an octant. The apparent high multipole excess is a 
consequence partly of the nonlinear corrections to the theory 
slightly underpredicting the nonlinear power, and partly of the 
shot noise due to the discrete nature of the particles. 

We computed the power spectrum of the lensed microwave 
background map and compared it with theoretical expecta- 
tions. Fig. [3] shows the excellent agreement between the the- 
oretical power spectrum from CAMB and the power spec- 
trum generated from the simulated lensed map. For complete- 
ness, we also show the power spectrum of the unlensed input 
map and the theoretical power spectrum used to generate it. 
The simulation confirms the theoretical prediction that lensing 
smoothes out acoustic features in the microwave background 
and moves power from large to small scales. Parenthetically, 
we would like to remark that the theoretical power spectra 
in Figs. 12] & [3] were generated with CAMB using the param- 
eter k_eta_max_scalar set to 20 x l_max_scalar. 
For k:_eta_max_scalar = 2x l_max_scalar (which 
is normally used) the n power spectrum, as well as the lensed 
microwave background power spectrum, become inaccurate 
at high multipoles; the n power spectrum is suppressed be- 
yond £ ^ 1000, and the lensed microwave background power 
spectrum beyond t ^ 4000. 

A particularly interesting quantity is the difference be- 
tween the lensed and the unlensed microwave background 
maps. Fig. |4] shows this quantity on a 7 degree x 9 de- 
gree patch projected onto a cylindrical-equal-area grid from 
the HEALPix sphere. The figure shows the striking large- 
scale non-Gaussian features that lensing induces on the mi- 
crowave sky. Lensing reconstruction techniques use these 
large-angular-scale correlations of small-scale features in the 
microwave background to recover the lensing potential. 

2.4. SZ Signal 

The SZ effect is a secondary distortion of the microwave 
background that results when cosmic microwave radiation 
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scatters against free electrons. It is commonly considered 
to ha ve two main components dSunvaev & Zeldovichlll970i 
[T972I) . The thermal SZ (TSZ) term arises from inverse 
Compton scattering of the CMB with hot electrons, predomi- 
nantly associated with shockheated gas in galaxy clusters and 
groups. The kinetic SZ (KSZ) term is a Doppler term com- 
ing from scattering with electrons having fast peculiar mo- 
tions. Another common distinction is that the KSZ effect has 
a nonlinear component associated with the small-scale intr- 
acluster medium and a more linear component coming from 
the large-scale intergalactic medium (IGM). The signal aris- 
in g from the latter was first calculated for the linear regime 
bv lOstriker & Vishniad ( [19861) and lVishniad (fT987) and is of- 
ten referred to as the OV effect. Below, we first write down 
the formalism for the nonrelativistic limit and then the more 
general relativistic case. 

Modeling the SZ effect requires knowing the number den- 
sity tie, temperature Te, and velocity of the electron distri- 
bution. In the nonrelativistic limit, the change in the CMB 
temperature at frequency v in the direction h on the sky is 
given by 

-^(n)=i-—) =f.y-b, (1) 

^CMB V ^CMB / tsz V ^CMB / ksz 

where the dimensionless Compton y and Doppler b parame- 
ters. 



y = j rieTedl = / ^^dr, 



b=^ I n^viosd/ = / Aosdr, 



(2) 
(3) 



are proportional to integrals of the electron pressure and mo- 
mentum along the line of sight (los), respectively. The dimen- 
sionless temperature, line-of-sight peculiar velocity, and the 
optical depth through a path length dl are given by 



kuTe 



= 1.96X 10" 



VkeV J 



Vlos 



Aos^ — = 3.34xl0-\ 

c VlOOkm/s 



Vlos 



dr EE ajriedl = 2.05 x 10" 



10"3 cm^ 



dl 
Mpc 



(4) 
(5) 
(6) 



respectively. For a typical cluster, the Compton y parameter 
is expected to be approximately an order of magnitude larger 
than the Doppler b parameter 

The nonrelativistic TSZ component has a frequency depen- 
dence as specified by the function. 



= xcoth(ji:/2)-4, x = /iz^/(^b7cmb), 



(7) 



which has a null at « 218 GHz. The distortion appears as 
a temperature decrement at lower frequencies and as an in- 
crement at higher frequencies relative to the null. The non- 
relativistic KSZ component is independent of frequency, but 
the sign of the distortion depends on the sign of the line-of- 
sight velocity. We chose the convention where vios > if the 
electrons are moving away from the observer. 
In the general, relativistic case, the change in the CMB tem- 
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Fig. 5. — Power spectra of the Sunyaev-Zel'dovich signals at 148 
GHz. Shown are the power spectra of the thermal SZ component (red), 
the kinetic SZ component (blue), and the full SZ signal (black). The 
soUd lines represent the model in this simulation, and the dashed hnes 
represent a mode l with no f eedba ck or star formation (adiabatic model 
IBode et a l. 2009; Trac et alj |2Q09|) . The green line shows the com- 
parison with the thermal SZ power spectrum derived analytically from 
Komatsu & Seljak (2002) and taken from the WMAP LAMBDA website 
( http://lambda.gsfc.nasa.gov/product/map/current/ 1. 



perature at frequency v in the direction n is given by 



Ar 

TCMB 



(«) = 



-pu^+OeCi+e^eCi) 



dr, 



(8) 



where t he F's and C's are known frequency-dependent coeffi- 
cients (INozawa et al.l il998). Note that Eq.[T]contains only the 
first-order terms in Eq. [8] and that fi^ = Yq. We use Eq. |8]to 
calculate the full SZ signal in these simulations. 

Sky maps of the SZ effect are made by tracing through 
the simulated electron distribution and projecting the accu- 
mulated temperature fluctuations onto a HEALPix grid with 
Nside = 8192. Using the saved simulation data described in 
^2.21 SZ maps are constructed with contributions from 3 com- 
ponents: massive halos at z < 3, low-mass halos and the inter- 
galactic medium at z < 3, and the high-redshift Universe for 
3 <z< 10. 

2.4.1. Gas model for massive halos at z <3 

The hot gas distribution associated with massive dark mat- 
ter halo s is cal culated using a hydrostatic equilibrium model 
(Osti-ikefeFaD^OOS; Bode et al. 2007, 2009). He re we sum- 
marize the general method and refer the reader to lBode et alj 
( 120091) for additional details, in particular on calibrating star 
formation and feedback against recent optical and X-ray ob- 
servations. In the large-scale-structure simulation, massive 
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dark matter halos with Mpop > 2 x 10'-' IC^Mq and z < 3 are 
postprocessed to include intrahalo gas and stars. These halos 
are expected to be the dominant contribution to the TSZ terms 
in Eqs.[I]and[8]for the scales of interest (1000 < / < 10000). 
Lower mass halos are accounted for differently, as we explain 
later. 

For each massive halo, the collection of N-body particles 
tagged by the FoF finder are enclosed in a nonperiodic, Carte- 
sian mesh that is centered on the center of mass, has one axis 
parallel to the line of sight, and a cell-spacing that is twice the 
N-body spline softening length. The particles are assigned 
using a cloud-in-cell scheme to define the matter density p,„ 
on the mesh. The gravitational potential (/) is then calculated 
using a nonperiodic Fast Fourier Transform (FFT), and the lo- 
cation of the potential minimum is chosen as the halo center 
This procedure preserves the concentration, substructure, and 
triaxiality of the dark matter halo. 

Assuming hydrostatic equilibrium and a polytropic equa- 
tion of state, the corresponding gas density pg and pressure Pg 
on the same mesh are given by 



Pg = PgP^P"'^-'\ 



r/(r-i) 



(9) 
(10) 



respectively, where the dimensionless temperature variable. 



i; = l- 



r-1 



^K'/'-^o), 



(11) 



is related to the gravitational potential. iBode et all (l2009h 
adopted a polytropic index F = 1 .2 (consistent with the value 
of F adopted by other workers and the value found in hydro- 
dynamic simulations) and determined the other two free pa- 
rameters for the polytropic model, namely the central gas den- 
sity pgfl and pressure Py o, by matching the gas mass fraction 
and temperature to recent X-ray observations ( Vi khlinin et al.l 
12006; Sun et al.,,2009). The electrons are assumed to have the 
same distribution of temperatures and velocities as the fully 
ionized gas. 

For calculating the KSZ terms in Eq. [T] and [8] the gas ve- 
locity on the mesh can be modeled in two ways. The veloc- 
ity field can be taken from the dark matter particles; this ap- 
proach allows for velocity substructure, although the disper- 
sion would be inconsistent with the hydrostatic equilibrium 
assumption. Alternatively, the entire gas distribution can be 
assigned a single, mass-weighted average peculiar velocity; 
this is the approach we adopt for self-consistency. Comparing 
the two options, we underestimate the temperature variance 
< AFj.^^^ > from massive halos alone by only ^ 4% when ig- 
noring velocity substructure. Correspondingly, the KSZ angu- 
lar power spectrum from massive halos alone is also underes- 
timated, but only by < 1% at / = 1000 and - 6% at / = 10000. 

When a given halo is added to the SZ maps, mass and mo- 
mentum are conserved by requiring that the baryonic mass 
be equal to the cosmic baryon fraction fib / fim times the halo 
mass MpoF- This is accomplished by first sorting the mesh 
cells by their gravitational potential, from most negative to 
least negative. The mesh is then traversed in order and cells 
are tagged until the accumulated baryonic mass, in gas and 
stars, reaches the mass limit. The outer radius of the gas is 
usually beyond the cluster virial radius. Remaining cells do 
not contribute any SZ signal. This procedure results in a tri- 
axial halo gas distribution with a physically-motivated cutoff. 

We find that the gas associated with these massive halos 
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Fig. 6. — The full Sunyaev-Zel'dovich signal at 148 GHz, including ther- 
mal and kinetic SZ contributions from both the galaxy clusters and the inter- 
galactic medium, in a ^ 10 degree X 10 degree patch. 

makes up approximately half of the TSZ flux, but it is the 
dominant contribution to the TSZ angular power spectrum. 
It accounts for ^ 90% of the total TSZ power for the scales 
of interest (1000 < I < 10000). On the contrary, this gas ac- 
counts for only ~ 15% of the total KSZ power because of the 
relative importance of the IGM to this signal, as discussed in 
j2A2l and E43] 

2.4.2. Low-mass halos and intergalactic medium at z <3 

Low-mass halos and the intergalactic medium at z < 3 are 
modeled using all other saved N-body particles not associ- 
ated with the massive dark matter halos discussed in ^2.4. II 
While we are able to identify halos down to Mpop = 6.82 x 
lO'^ /i~'M0 (100 particles), there is insufficient resolution 
to accurately apply the hydrostatic equilibrium model. In- 
stead, an alternative approach motivated by the virial theorem 
is adopted, one that can also be extended to particles not as- 
sociated with identified halos. The intergalactic medium is 
expected to contribute only a small fraction of the TSZ signal, 
but it is very important to include it for calculating the KSZ 
signal. 

According to the virial theorem, the internal energy of a 
virialized halo is twice its gravitational potential energy. Dur- 
ing the formation of a halo, the dark matter undergoes violent 
relaxation and the gas gets shock-heated. The gas tempera- 
ture can be estimated directly from the halo mass and used 
to model the TSZ signal (e.g. .Cohn & Whitell200 9). but we 
chose to calculate the gas temperature from the velocity dis- 
persion of the dark matter. This approach is then applicable 
to all particles rather than just ones in identified halos. 

For each particle, the local density p„, and velocity disper- 
sion (T,, are calculated using 32 neighboring particles. For cal- 
culating the SZ terms, the effective pressure and momentum 
of the gas is taken as 



rim 



(7- l)p„,a^ 



(12) 
(13) 
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where 7 = 5/3 is the ideal gas adiabatic index and /i = 4/(3 + 
5Xii) = 0.588 is the mean atomic weight for a fully ionized 
gas. Helium reionization is assumed to occur instantaneously 
at z = 3. A temperature floor of lO"* K is imposed for the 
photoionized intergalactic medium. 

When tracing through a given particle and projecting onto 
a HEALPix grid, the effective solid angle subtended by the 
particle is taken into account. A nearby, low-density particle 
can span a number of map pixels and simply assigning it to 
a single pixel will result in a large and incorrect temperature 
fluctuation. Mass and momentum are also conserved. Since 
the entire gas distribution at z < 3 is accounted for, the maps 
have a smooth transition from the concentrated halos to the 
diffuse intergalactic medium. 

We find that the low-mass halos and the IGM make up ap- 
proximately half of the TSZ flux, but contribute to the angular 
power spectrum only by ^ 10% at I = 1000 and ^ 15% at / = 
10000. This is in agreement with Hallman et al.l (12007 ) who 
used a high-resolution hydrodynamic simulation and found 
that while roughly one-third of the TSZ flux comes from low- 
mass halos (M < 5 X IO^^Mq) and the IGM, these compo- 
nents add_<_10%to the TSZ angular power spectrum (for erg = 
0.9). iHallman et al.l ([2009) also found that the angular power 
spectru m of the filamentary warm-hot intergalac tic medium 
(WHIM ICen & Ostrike^ [19991: iDav e eTaP '200B is only < 
1% of the total. Furthermore, iHer nandez-Monte agudo et aP 
([2006") found that while lower density gas (5 < 10) makes up 
approximately half of the baryon budget, it gives rise to only 

5% of the TSZ flux in the local universe. 

For low-mass halos in the mass range 6.82 x 10'^ < 
MpoF {h'^Mo) < 2 X 10'^, there are adequate particle numbers 
(100 to 282) to estimate the halo velocity dispersion. The TSZ 
signal from this component is sufficiently well approximated 
with the above method, considering that these halos contribute 
< 10% to the TSZ angular power spectrum. However, for fila- 
ments with typical overdensities of 10, there are only 10 parti- 
cles per comoving (/!"'Mpc)^. While it is difficult to calculate 
the velocity dispersion accurately, an approximate solution is 
acceptable since the filaments contribute only < 1 % to the an- 
gular power spectrum. 

The KSZ signal from the low-mass halos and the IGM are 
relatively more important. We find that at / = 1000, they ac- 
count for ^ 15% of the total power, comparable to the con- 
tribution from massive halos, and at / = 10000, they account 
for ^ 40%. The majority of the KSZ signal is fro m the h igh- 
redshift universe (3 < z < 10), which we discuss in 32.4.31 The 
KSZ signal is relatively less sensitive to resolution since the 
bulk peculiar momentum is well captured in the large-scale- 
structure simulation. While velocity substructure is not well 
resolved, we have already demonstrated in 32.4.1l that ignor- 
ing it only results in percent level differences. However, given 
the modest resolution of the large-scale-structure simulation, 
we caution against using the maps for detailed studies of the 
SZ signal from the WHIM and IGM. 

2.4.3. High-redshift Universe at 3 < z < 10 

The high-redshift universe at 3 < z < 10 is modeled us- 
ing the redshift shells containing the surface density fields for 
mass S,„ and line-of-sight momentum S,„,,. At these higher 
redshifts, it is very reasonable to assume that the gas traces the 
matter distribution because the scales of interest are mostly 
still in the linear regime. Furthermore, since both the nonlin- 
ear mass and the collapsed mass fraction above a fixed mass 



rapidly decrease towards higher redshifts, the typical temper- 
ature of the gas is predominantly set by photoionization rather 
than shockheating. 

For each redshift shell, the electron temperature, line-of- 
sight peculiar velocity, and column density are taken to be 



(14) 
(15) 



Te = 


10'*K, 






Vlos ■■ 








(a^) ( 









(16) 



Hydrogen reionization is assumed to occur instantaneously at 
z = 10 and helium is only singly ionized. The Thomson opti- 
cal depth for electron scatteri ng is tt ~ 0.08, con sistent with 
the WMAP 5-yeai- results dlDunklev et alj 120091) . Since the 
electron distribution is warm rather than hot, the SZ temper- 
ature fluctuations are nonrelativistic and Eq.[8]reduces to Eq. 
[U There is only a very small contribution to the TSZ signal 
from the warm electron distribution. The average Compton y 
parameter is only ^ 10"^ over this redshift range. However, 
there is a substantial contribution to the KSZ signal because 
of the large, coherent velocities in the IGM. 

Since the light cone is constructed by replicating and stack- 
ing from the periodic simulation box without random rota- 
tions, this procedure introduces excess power in the KSZ sig- 
nal at / < 500. The excess power shows up in the intergalactic 
medium component, but not in the halo component previously 
discussed in 32.4.11 Therefore, the intergalactic medium com- 
ponent of the KSZ, from the full redshift range < z < 10, is 
filtered before adding it to the final SZ maps. A simple filter 
w(l)= l-exp[-(Z/500)^] is chosen to suppress the large-scale 
excess, but preserve power at / > 1000. The remaining low-Z 
power is now more consistent with the level expected for the 
large-scale OV signal. 

Since the simple filtering modifies the signal at I < 1000, 
the maps should not be used to predict the KSZ signal at these 
scales. Furthermore, we have neglected to model the con- 
tribution from the inhomogeneous reionization epoch. This 
component is expected to be larger than the contributio n from 
the reionized universe (z 5, 6) at I < 1000 (e .g. , Santos et al.l 
[2003t iMcOuinn et al.l[200l iZahn et al.l[2005[) . The KSZ sig- 
nal on these larger scales is subdominant compared to the 
lensed microwave background and to the TSZ signal at fre- 
quencies away from the null. 

2.5. Infrared Point Sources 

Star formation in the Universe leads to the production 
of dust grains which absorb the ultraviolet radiation in 
the environment and re-radiate it in the infrared range. In 
the context of microwave background observations, the 
redshifted infrared emission of these sources constitutes an 
important contaminant on small angular scales at frequencies 
higher than ^150 GHz. In this section, we attempt to model 
the angular anisotropics that infrared (IR) galaxies introduce 
in the microwave sky by combining the cosmological simu- 
lation presented above with a model for the infrared galax y 
population, which is partially based upon lRighi et al.l (120081) . 



2.5.1. The basics of the infrared model 
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We first approximate the spectral energy distribution (SED) 
of the infrared galaxies by a modified black body of the form 

oc (j)(ujaust) = [Biy(7dust)-Bi/(7Lmb)] , (17) 

where /3 is a spectral index, Tcmb is the primary microwave 
background monopole temperatu re, and Tjust i s the dust tem- 
perature evolving in redshift as in lBlainl([l999h : 



7dust(2) = 



7;r+ca(i+^)'^^-i) 



l/(4+/3) 



(18) 



Here Tq is the dust temperature if the galaxy was at z = 0, and 
the symbol B^{T) denotes the black body spectral intensity at 
temperature T. As we shall see below, the free model parame- 
ters (3 and Tq will determine the spectral index for the effective 
power law describing the change of the infrared intensity ver- 
sus frequency, I^, oc i^" (i.e., a = a(/3, Tq)). 

The star formation rate (SFR) in a given galaxy seems to 
be line arly correlated w ith the dust luminosity (Kennicutt re- 
lation, I^nnicut3[T998'). At the same time, there is observa- 
tional evidence that the SFR reached a maxim um at epochs 
corresponding to z ^ 2-4 jHopkins & Beacomli2006l) . These 
two factors combined suggest introducing a redshift window 
function for the infrared luminosities in galaxies close to the 
measured SFR behavior. In our model, we adopt the depen- 
dence 



Ws, 



,(z) cx exp(-l/z) exp (- (z/3.5)') . (19) 



This window however gives more weight to the re dsfhift 
range z G [ 1 , 2] if compared to the data of Hopkins & B eacomI 
([200^, due mainly to the absence of halos in our simulation 
atz> 3. 

In practice, given the simulation halo catalog, our approach 
consists of populating each halo with a number of galaxies 
with a given infrared luminosity. We adopt a mass range lim- 
ited by Mi,M2 (Ml < M2), in such a way that halos below 
Ml will host no infrared galaxies. Halos with masses M E 
[Mi,M2] will typically host one infrared galaxy of total in- 
frared luminosity L/r(M, z) = Wsfr(z) ■ W^ooiiM) ■ U ■ (M/M2). 
For more massive halos, the average number of infrared 
galaxies will be given by M/M2, each of them with a total in- 
frared luminosity given by Lir(M,z) = WsFRiz) ■ Wcooi(M) ■ L^,. 
The function WmoiiM) = exp{-M /Mc„oi) accounts for the fact 
that bigger halos need longer times to cool down before being 
able to form stars. L^, defines a luminosity parameter The 
parameters Mi, M2, and Mcoi, combined with the spec- 
tral parameters /? and Tjust introduced above, will define our 
infrared model. In the simulations, the actual number of in- 
frared galaxies present in a given halo is driven by a Poisso- 
nian realization of its average number. We also introduced 
an effective scatter of 15% in the mean adopted value for 
the dust spectral index and in the infrared luminosity of each 
galaxy, partially motivated by the scatter in the effecti ve spec- 
tral index of the intensity found in the analyses bv ,Knox et al.l 
(l2004 . 

2.5.2. The angular power spectrum 

Below we calculate the angular power spectrum of our in- 
frared model theoretically. This allows us to approximate the 
power spectrum efficiently, after adjustments to the model pa- 
rameters, without generating HEALPix maps. For each IR 
galaxy, the spectral luminosity can be computed from the 



bolometric luminosity and the adopted form of the SED: 
Lir{M,z) 



J dv cj)(i', Tdusi) 



(20) 



where Tdust is a function of redshift (Eq. [TsT l. For a flat Uni- 
verse, the measured spectral flux is therefore given by 



" 47rr2(l+z)' 



(21) 



where v = i^'/(l + z) and r is the comoving distance to the IR 
galaxy placed at redshift zif). The average spectral intensity 
generated by the IR galaxy population can then be written as 
the sum of the contribution of all galaxies per unit solid angle: 



T ; dn a(f)Lu'{M,r) 
I.= Idrr' I dM—ir)N,iM)^2^^ 



(22) 



The comoving halo mass function is given by dn/dM, and 
Ng{M) provides the average number of infrared galaxies in a 
halo of mass M. The cosmological scale factor is expressed 
by a(r), and the typical spectral luminosity for each of the 
infrared galaxies in the halo is given by Ly/(M, r), with i^' 
being the frequency observed at reds hift zj r). Li,i{M, r) is de- 
termined by the model described in ^2.5. H and Eq. |20l 

However, we are interested in the angular fluctuations of 
the infrared spectral intensity, and this is linked to the angu- 
lar variation of the number of infrared sources. The number of 
infrared galaxies is driven a priori by Poisson statistics. How- 
ever, halos tend to cluster in overdense regions of the large 
scale matter density field, and this introduces another modu- 
lation in the infrared galaxy number density that is known as 
the correlation term (iHaiman & Knoxl2000l : ISong et alJ2003l : 
Righi et al. 2008)). In terms of the spectral intensity 1^, it is 
easy to show that the Poisson term for the angular anisotropics 
enerated by the i nfrared galaxy population can be written as 



generated by 
ilRighi et all 



mm 



drr' 



dn ( a(r)L,^'(M,r) 

dM —(r)N,{M) '/\ 
dM ^ \ A-Kr^ 



(23) 



whereas the correlation or clustering term depends upon the 
halo power spectrum: 



%Pm(k,r) 



dM^{r)b{M,r)Ng{M) 



a{r)L^>(M,r) 



n 2 



47r 



(24) 

where ^ = (/ + 1 /2)/r. The bias factor b{M, r) expresses how 
halos are distributed compared to dark matter, in such a way 
that b^(M, r) is defined as the ratio of the halo power spec- 
trum and the matter po wer spectrum (Pmik, r)). Here we have 
adopted the model by ISheth & TormenI ( 1199 9^. which pro- 
vides a good fit (few percent accuracy) to the simulated halo 
power spectrum. The last two equations are obtained after us- 
ing the Limber approximation , which sho uld be accurate to 
better than 2% for / > 20 (e.g.. lVerde et"ani2000l) . 

2.5.3. Observational constraints 

The parameter space is swept in such a way that our infrared 
model satisfies the following observational constraints: 

• T he Cosmic In f rared Background (CIB), as estimated 
in iFixsen et alj ( 119981) from COBE/FIRAS, is the up- 
per limit for all total infrared intensity that our model 
provides. 
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Fig. 7. — Left panel: Cumulative source counts of our infrared model at 353 GHz versus flux (red triangles) . The solid line corresponds to the fit to the SCUBA 
data at the same frequency iCoppin et al. 2006). The high flux sources match the observed BLAST counts (Patanchon et al. 2009 and BLAST team, private 
communication). Middle panel: Contribution of our infrared source population to the total cosmic infrared background observed by COBE/FIRAS I Fixsen et al.l 
[1998). We show the three observing frequencies of ACT plus 353 GHz. Dashed lines indicate 68% confidence regions. Right panel: Theoretical estimate of the 
angular power spectrum of our mock infrared source population. The thick black solid line con'esponds to the clustering term , and the thin line to the Poisson 
term. The total is given by the red line, which is below the upper limit suggested by ACBAR at / = 2600 iReichardt et alJ2009l) . 



• The choice for spectral parameters (3 and Tdust must 
be such that the effective spectral index for the spec- 
tral intensity {Ii, oc z^") is close to a w 2.6 in th e range 
V G [145,350] GHz. This value is found in .Knox et alj 
(2004), and is bas ed on the projection to high redshifts 
of the analyses of lDunne et alJ (120001) of a low redshift 
infrared galaxy sample. Nearby galaxies show harder 
spectral indexes, while high z sources are closer to the 
peak of the dust SED and hence provide a shallower 
slope. The quoted value of the spectral index is prac- 
tically equivalent to the slope of the CIB at the same 
frequency range. Thus, imposing a « 2.6 makes the 
infrared background generated by our infrared sources 
have a slope close to that of the CIB in the frequency 
range i/G [145,350] GHz. 

• The source counts of our infrared model must be 
compatible to the s ource counts provided by SCUBA 
(ICoppin et al.ll2006l) in the flux range S^s^ghz G [1j50] 
mJy. These source counts are in agreement with 
a number of other su bmillimeter observations (e.g. 
lAustermann et al.ll2009h . 

• The source counts of our infrared model must be com- 
patible with the high flux in frared source population 
recently seen by BLAST ( Patanchon et al.l l2009l and 
BLAST team, private communication). 

• The clustered component of the the angular power spec- 
trum should dominate over t he Poisson term for / < 
1000 as observed by Spitzer dLagache et al.ll2007l) and 
BLAST (Vieroetal. 206§. 

• The angular power spectrum introduced by our IR 
source population should not exceed the upper limit 
derived from ACBAR observations at 150 GHz after 
fixing to the WMAP 5-year best-fit value with a 
running spectral index (erg = 0.81), i.e. I^Ci/{2tt) < 
3liiiKf at / = 2600 dReichardt et al.ll2009l) . 



2.5.4. Final model 

We modeled two different source populations: the basic 
population essentially defines the properties of the whole IR 
model, whereas the very high flux population matches the high 
flux IR source counts detected by BLAST (BLAST team, pri- 
vate communication). Each of these two models requires a 
different choice for the parameter set, as shown in Table [T] 

In order to avoid the limitation imposed by the minimum 
mass of the halo resolved in our simulation, we added a 
mock catalog of halos of masses between lO** Mq and the 
minimum resolved halo mass, Mm,„ = 6.82 x 10'^ /i~'Mo = 
9.61 X 10'^ Mq. In order to do this, we assumed that the 
mass functio n in th is mass range was dictated by the fit of 
iJenkins et al.l (120011) . We used the population of low mass 
halos resolved in the simulation as tracers of the mock ha- 
los appended in this low mass range. We divided the mass 
range [lO** Mq,M„„„] in logarithmically spaced bins, and, in 
each of these bins, computed the average number of halos per 
resolved halo that is found in a thin mass interval centered 
at Mc = 1.10 X lO'-' Mq. That is, for each resolved halo of 
mass close to M^, we computed the average number of ha- 
los present in each of those low mass bins, and subsequently 
made a Poissonian realization of this average. The masses 
of these new mock halos are randomly distributed uniformly 
within the corresponding mass bin width. These mock halos 
were randomly distributed (with a Gaussian profile) in spheres 
of roughly 100 times the parent (resolved) halo virial radius, 
which should correspond to ^ 10 Mpc. Therefore, for each 
real halo of mass close to M^ we populated a cloud of lower 
mass halos in an sphere of roughly 10 Mpc radius and cen- 
tered on the parent halo itself. This procedure provides our 
simulation with a low mass halo catalog that is not Poisson 
distributed, but rather is correlated to the population of small 
but resolved halos close to M = Mc. In practice, only mock 
halos of mass above Mi were generated, since they are the 
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TABLE 1 

Parameters for the two populations considered for the infrared 

POINT SOURCE model 



Model 


/9 


Tq 


Ml 


M2 


Mcool 








(K) 


(Mq) 




(Mq) 




Basic 


1.4 


25 


2.5 X 10" 


3 X lO'^"* 


5 X 10"* 


1.25 X 10*2 


Very High Flux 


1.3 


40 


5 X lO"* 


1 X 10'5 


8 X lO"* 


6x 10"* 



ones eventually hosting IR galaxies. 

The very high flux population consists of few 300 in the 
whole octant) infrared galaxies above 200 mJy at 353 GHz. 
This very high flux population models the tail observed by 
BLAST at 250, 350 and 500 /im (BLAST team, private com- 
munication). The high flux tail should be at or above the 10 
mJy level in the 145 - 220 GHz channels, contributing sig- 
nificantly to the power spectrum. Some of this very high flux 
population includes sources believed to be lensed and magni- 
fied by the intervening structure (Lima et al. 2009). 

The total cumulative source counts above a given flux 
threshold for the sum of the basic and the very high flux 
populations are displayed in the left panel of Fig. |7] We 
show fluxes above 1 mJy (which is about SCUBA's detection 
threshold), although our IR source population extends to dim- 
mer fluxes. It is actually the numerous but dim IR source 
population that is responsible for most of the CIB. In our 
model, our IR source population generates between 50% and 
80% of the total CIB detected by FIRAS (middle panel). This 
is consistent wit h BLAST observations ( iMarsden et al.ll2009l : 
iDevUn et alj|200 9), which indicate that most, and perhaps all, 
of the CIB is generated by infrared sources. 

In the right panel, we display the expected angular power 
spectrum of our total IR source population. At very high 
multipoles the Poisson term must dominate, whereas the clus- 
tering term takes over at / ^ 1,000. Again due to the con- 
struction of the light cone by replicating and stacking the pe- 
riodic N-body simulation box without random rotations, ex- 
cess power is introduced in the infrared point source signal 
at / < 300, as was the c ase fo r the kSZ signal from the inter- 
galactic medium (see § l2.4.3l l. We do not filter this signal to 
suppress the low-/ excess as we did for the kSZ, as this signal 
consists of discrete sources instead of a smooth component, 
and such a filter would disturb the real space point source dis- 
tribution. For / > 300, the theoretical expectation of the angu- 
lar power spectrum coincides with the actual power spectrum 
from the simulated infrared map. 

The particular values for the model parameters for each of 
the two populations are displayed in Table [T] Note that the 
masses in this table refer to friends-of-friends halo masses. 
For the basic population, the values of the dust tempera- 
ture and infrared luminosity are not far from those found by 
iDve et al.1 (12009) when comparing BLAST infrared sources 
to radio and mid-infrared counterparts in the Chandra Deep 
Field South. For the few galaxies in the very high flux pop- 
ulation, characteristic masses and infrared luminosities must 
be increased by roughly an order of magnitude, whereas the 
typical dust temperature increases from 25 to 40 K. 

2.6. Radio Point Sources 

Our model for radio point sources foUows in spirit the treat- 
ment of the subject by IWilman et al.l (120081 hereafter W08). 
Similar to W08,we separate the radio-loud active galactic 
nuclei (AGNs) into two populations, a high-power family 



that evolves strongly with redshift, and a low-power one 
with mild redshift evolution. These two populations can be 
roughly ide ntified with the type II a nd I classification of radio 
sources by Fanaroff & Rilev (1974), respectively (however, 
see discussions in Willottet al. 2001). We also follow W08 
and iJackson & W all (1999) to employ a relativistic beaming 
scheme to separate the contributions to the total flux from the 
compact core and the extended lobes of radio sources. 
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Fig. 8.— The 151 MHzR LF at z < 0.3 constructed from the 3CRR survey, 
using data from lLaing etan fl983). The blue squares are for FR lis, and the 
green triangles are for FR Is. The red pentagons are the model prediction 
(sum of the two populations). 

Unlike W08, we do not consider radio-quiet AGNs and 
star-forming galaxies in the radio source model. As our pri- 
mary goal is to simulate the microwave sky at relatively bright 
flux limits (e.g., > 0.1 mJy), ignoring these populations has 
negligible effects. Furthermore, the infrared source model- 
ing presented in ^2.5l is a self-contained model for the star- 
forming galaxies, and therefore the combined infrared and ra- 
dio source models should adequately account for the majority 
of sources of interest to the current generation of microwave 
experiments. 

Our model is a way to populate dark matter halos with ra- 
dio sources in a Monte Carlo fashion, according to simple 
prescriptions on (1) the spectral energy distribution of the 
sources, (2) the number of radio sources as a function of halo 
mass and redshift, and (3) the (radio) luminosity distribution 
(LD) of the sources. Below we outline the most important 
aspects of the model. 
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2.6. 1 . Compact versus extended emission of radio sources 

The extended lobes of radio sources are characterized by 
a steep falling spectrum. On the other hand, the emission 
from the core region of a radio source is usually flatter and 
can ex hibit a more complex spectral shape (see e.g. jLin et alj 
l2009bi) . In addition, when the jet axis is close to the observer's 
line-of-sight, the fluxes from the core (and/or the base of the 
jets) will be Doppler-boosted. 

Following W08, we assume the observed core-to-lobe flux 
ratio, Rohs, is related to the intrinsic value, Ri^u via ^obs = 
RintB(9), where 9 is the angle between the jet axis and the 
line-of-sight, 8(6) = [(I- (3 cos 6)-^ + (I + ^ cos 9)-^]/ 2, 13 = 
7~^, and 7 is the Lorentz factor of the jet. 

Given the intrinsic luminosity Lj^ = i-c.int+A.int of a source, 
where the subscripts c and 1 refer to the core and lobes, the 
beamed luminosity is then 



Lb, 



)eam 



— ^c,beam "t" ^l,int 
= (1 +^obs)i'l,int 
= LUl+ Robs) 1(1+ Rind 



(25) 



We use empirically derived values of /?int and 7 from W08, 
and we assume a uniform distribution of cos 9 in order to sep- 
arate Lcbeam and Lijnt givCn Lbeam- 

Since the compact core emis sion is usually observed to ex- 
hibit a complex spectral shape (iLin et al]|2009b<) . we model it 
with a third-order polynomial, log/^ = "Y^. a,(logz^)', where v 
is in GHz. We assume the extended sourc es obey fy oc v~^- ^, 
which is also consistent with observations (iLin et al.ll2009bl) . 

2.6.2. Halo occupation number and the 151 MHz RLF 

Our first task is to reproduce the local 151 MHz radio 
luminosity function (RLF), which is effectively a convolution 
of the halo mass function, the halo occupation number, 
and the radio luminosity distribution. The main reason 
for targeting the RLF at 151 MHz is because at such low 
frequencies, one can safely assume that the fluxes are mostly 
dominated by the steep spectrum, extended sources, and the 
resulting RLF is not very sensitive to biases due to orientation 
effects. In addition, this RLF is not greatly affected by 
synchrotron/inverse Compton losses ( Blundell et aLl ri999l) . 
For both FR 1 and FR 11 populations, we assume their halo 
occupation number (HON) is given by N(M) = Nq(M/Mq)°'. 
As for the LD, we assume it is independent of halo mass, 
and adopt two-segment power-laws that roughly mimic the 
observed RLF: 



p(L): 



(L/LiT if L> Lb 
(L/Lh)" otherwise 



(26) 



We explore the parameter space spanned by the HON and 
LD for both radio source populations, until the combined RLF 
fits the observed one at z ~ 0.1, which is based on the 3CRR 
survey (Laing et al. 1983). A comparison of the model pre- 
diction and the observed 151 MHz RLF is shown in Fig. [8] 

2.6.3. Redshift evolution and source number counts 

Once our model successfully reproduces the low redshift 
151 MHz RLF, the radio sources are divided into two popula- 
ti ons, roughly correspondi ng to the type 1 and II classification 
of lFanaroff & RilevI (11974 ). These populations are given dif- 
ferent number density evolutions with redshift, so that both 




S (Jy) 

Fig. 9. — Source counts at four frequencies (1.4, 33, 95, and 145 GHz). The 
model prediction (one realization of a,'s - see Table|2) is shown as red solid 
points. All other symbols are observed source counts. The only exception 
is the curve in the 95 GHz panel, which is itself a model from de Zotti et aQ 
120051) , as shown inFSadler et al. (200 8). The da t a us ed in the 1.4 GHz panel 
are from Katgert^t alj a988i) : iWindhorst et alj fT993): Bondi et al. l2003p; 
Huynh et al. 1 2005). The data for 33 GHz is from Wright et al. 1 2009|; 
Clearv et alj I l2005t) . The data for 145 GHz is from ACBAR tReichardt et aP 
2OO9I1 

the high-z RLFs and source counts at 151 MHz can be fit- 
ted. For FR 1 sources, we consider an evolution of the form 
(l+z)^ up to Zp (and constant at z > Zp); for FR lis, we find 
that a much stronger redshift evolution is needed and, for con- 
venience, adopt an asymmetric Gaussian centered at Zp, with 
dispersions ct; and ct/, at z < Zp and z > Zp, respectively. (Note 
that there are many more FRls and FRlls at high redshift; 
for example, the comoving density of model FRll sources is 
^ 200 times higher at z^ = 1 .3 compared to at z 0.) In our 
model, Zp is a free parameter, which is adjusted to match ob- 
servations. 

Having matched the high-z 151 MHz RLF, we assume each 
source at 151 MHz is described by Lbeam = Lcbeam + Li.int, 
where the lobe component dominates. Given Ri„i and 7 from 
W08, and cos 6* picked from a uniform distribution, we de- 
termine Lc.beam- We ncxt adjust the shape of the core SEDs 
(by varying the coefficients a,) so that the source counts at 
f > 151 MHz are reproduced. These source counts at higher 
frequencies are dominated by the radio cores. In Fig. 9, we 
show a comparison of our model prediction (red solid points) 
with observed source counts at 1.4, 33, 95, and 145 GHz (see 
caption for observation references). In general the model pro- 
vides a good fit to the data'^. We note that for the 95 GHz 
panel, we are comparing our model to the predictions of the 
de Zotti et al. (2005) model. 

There are 24 free parameters in the model, which fit 24 ob- 
servational constraints (12 RLFs at various frequencies and 
redshifts, and 12 source count observations at 12 frequencies). 
(Note that qq, ai, 02, 03, and a are set to be the same for both 

We have also compared our source count predictions with observations 
at 0.15, 0.6, 2.7, 4.9, 8.4, 15, 20, 41, and 61 GHz, and found good agreement. 
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TABLE 2 

Model Parameters for Radio Sources 



Type 
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FRI 


j(,-2.6 
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(-0.12,0.07) 


(-0.34,0.99) 


(-0.75,-0.23) 
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4x 10'3 


0.1 


1024 


-1.55 





3 


0.8 


FRII 


j(,-2.8 
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(-0.12,0.07) 


(-0.34,0.99) 


(-0.75,-0.23) 


0.015 


3 X 10'^ 


0.1 


1027.5 


-1.6 


-0.65 




1.3 0.4 0.73 



^ We draw random urriform variates irr this range when assigning the spectral shape for cores. 
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Fig. 10. — The .simulated Galactic dust emission in a 10 degree X 10 
degree patch. 



FRJs and FRIIs.) Several of the parameters are degenerate; in 
Table |2] we present a parameter combination that can repro- 
duce the results presented here. 

For a more complete description of the model and fur- 
ther c omparison with available observational constraints see 
lLinetal.l(l2009ai) . 



3.1. General Map Properties 

Below we list some general properties of the maps them- 
selves. The full-sky HEALPix maps for this simulation 
are in the Ce lestial coordinate system with RING ordering 
riorski et al.H2QQ5h . The conversion between the HEALPix 
9 and (f) coordinates and right ascension and declination is: 
= ra, 90-6' = dec. We store all maps in units of flux per 
solid angle (as Jy/sr). To convert to temperature (as /xK) we 
use the conversion: 



dT _ 
dl Iksv^ 



(27) 



where x = hv/ksTcMB- To convert from SI to more convenient 
units requires /iK sr/Jy = 10^" K sr m^ Hz / W. When cal- 
culating the power spectrum with Polspice''' or anafast'"^, a 
correction for the pixel window function is required. 

The properties of the halos that are specified in the halo cat- 
alog are redshift, ra, dec, position of the halo potential min- 
imum, peculiar velocity, friends-of-friends mass, virial mass, 
virial gas mass, virial radius, integrated Compton y parameter, 
integrated kinetic SZ, integrated full SZ at the six frequencies 
(the latter three integrated within the virial radius, R200, and 
/?50o), the stellar mass, and the central gas density, tempera- 
ture, pressure, and potential. In the radio and infrared galaxy 
catalogs, the properties given for each source are redshift, ra, 
dec, and the flux in mJy at the six frequencies. 

The large-scale structure octant is mirrored into full-sky 
maps as follows: 



2.7. Galactic Emission 

For the contribution of Galactic t hermal dust emiss i on, w e 
use the "model 8" prediction from iFinkbeiner et alj ( 119991) . 
an extrapolation to m icrowave frequencies of dust maps from 
iSchlegel et all (11991 . We use the software included with 
these maps to interpolate onto an oversampled HEALPix grid 
with A^side = 8 192. This model is a two-component fit to IRAS, 
DIRBE, and FIRAS data, and is shown by Bennett et al] 
(j2003) to be a reasonable template for dust emission in the 
WMAP maps. 

The Galactic synchrotron and free-free emission is compa- 
rable to the Galactic dust emission at ^ 70 GHz and becomes 
dominate/ su bdominant at lower/higher frequencies (see e.g. 
ISmood[r998 ). Thus these signals are not expected to be a sig- 
nificant contaminant for ACT or SPT (which surveys target 
regions away from the Galaxy), although they will be impor- 
tant to consider for the Planck low-frequency channels. We 
leave the inclusion of these signals to lat er work, and note 
work by, for example, iLeach etaP (l2008h . which discusses 
the inclusion of these signals in microwave simulations for 
Planck. 



On = 00 (28) 

J. J. 

9n = <PQ + n- 

TT TT 

(Ps = 2~'Po + "2 

where n = 0, 1,2, and 3, and the subscripts and S refer to 
the northern and southern hemispheres of the map. The sub- 
script refers to the coordinates in the original octant. For 
the power spectrum analysis, this method for repetition in the 
Southern hemisphere is preferable to reflecting the North. Al- 
though simple mirroring eliminates real space discontinuities 
across the equator, it makes the simulated sky an even func- 
tion with respect to the polar axis. In the spherical harmonic 
transform, the associated Legendre functions are even or odd 
depending on /. Integrating mirrored hemispheres over them 
yields zero for odd-Z harmonic coefficients. Because rotating 
the sky transforms harmonics with common / into each other. 



3. DISCUSSION 



http://www.planck.fr/articlel41.html 
http://healpix.jpl.nasa.gov/html/facilitiesnode6.htm 
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odd-/ harmonics will be zero regardless of the orientation of 
the mirroring. 
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Fig. 1 1 . — Power spectra of the lensed microwave background (black), ra- 
dio galaxies (blue), infrared galaxies (green), and full Sunyaev-Zel'dovich 
signal (red) in this simulation at 148 GHz. The magenta line through 
the blue radio sour ce power spectrum gives the radio source power from 
the IToftbl atti^t^n ( 1998) model scaled by 0.64, which was found by 
IWright et al.l )2009fl to be a good fit to the WMAP 5-year source counts. 



3.2. Power Spectra of Simulated Components 

In Fig. [TT] we show the power spectra of the lensed mi- 
crowave background, infrared galaxies, radio galaxies, and 
the full SZ signal from the simulated maps at 148 GHz. The 
SZ power peaks at around / ~ 3000 and the clustering of the 
infrared galaxies is appar ent for / < 1000, co nsistent with 
measurements by Spitzer ( iLagache et al.ll2007h and BLAST 
(lyiero et al. 2009). No infrared or radio sources have been 
removed in creating the power spectra, which were made 
with Polspice on the full-sky maps. The power from the in- 
frared and radio sources, without masking out bright sources, 
is comparable to the lensed microwave background between 
I - 2000 and / ^ 3000. In Fig. [IT] the magenta line through 
the blue radio source power spectrum gives the radio source 
power suggested by the Toffolatti et al. (1998) model scaled 
by 0.64, which was found bv lWright et al.l (12009b to be a good 
fit to the WMAP 5 -year source counts. The two power spectra 
are in excellent agreement. 

3.3. Y-M Relation and Scatter 

We investigate the SZ flux versus cluster mass scaling re- 
lation for this simulation below. An accurate calibration of 
this relation off observations is important for assessing the se- 
lection function of a potential SZ selected cluster sample via 
simulations. The gas physics input into the simulations will 
determine the surface brightness of a given cluster based on 
its mass and redshift, which in turn determines how likely it 
is to be detected by a given observation strategy. 



To characterize the SZ flux - cluster mass scaling relation in 
this simulation we compare the integrated Compton y param- 
eter to cluster mass for the halos in the one unique octant. The 
integrated Compton y parameter is given hyY = (J ydil)d^{z), 
where y is defined in Eq. |2] dft is the solid angle of the clus- 
ter, and d^(z) is the cluster's angular diameter distance. For 
the self-similar model, we expect a virialized halo of mass 
Mvii to have a virial temperature of 



Tvi, cx [Mv„-£(^)]2/^ 
where for a flat ACDM cosmology 

E(z) = [n„(i+z)'+nA]'^^. 



(29) 



(30) 



If galaxy clusters were isothermal, then the Compton y pa- 
rameter would satisfy the relation 



F Oc/gasMhalor, 

where /gas is the cluster gas mass fraction, 
above equations, we would expect 



(31) 

Combining the 



(32) 



for the self-similar relation. Deviations from the self-similar 
relation are not unexpected as clusters are not isothermal, not 
always in virial equilibrium, and /gas need not be constant 
with variations in cluster mass and redshift. 

The gas physics in this simulation has been calibrated off 
of X-ray observations of clus ter g as f ractions as a function of 
temperature as discussed in ^2.41 and iBode et al] ([2009). For 
each halo with Mfof > 2 x 1 0'^^ Ii'^Mq and z < 3 (whose gas 
model is described in ^2.4.11 ) we calculate ¥200 and isoo, which 
are the projected Compton y parameter in a disk of radius /?200 
and /?50() respectively. Here /?2oo and /?50() are the radii within 
which the cluster mean mass density is 200 and 500 times 
the critical density at the cluster redshift. We then rank the 
clusters by mass and bin them in mass bins of Alog(M2()o) = 
0.05 and AIog(M5oo) = 0.1 for M200 and M500 respectively. 
The mean of Ma as well as the mean and standard deviation 
on the mean of Fa x E(z)~^^^ are then calculated for each bin. 
These values, forM2oo and Y200, are plotted in Fig. [12] along 
with all the individual halos with M200 > 2 x lO'^^M© in the 
octant. 

The binned values are then fit to the power-law relation 



Ya 



Eizp^ 



= 10^( 



Ma 
lOi^M^ 



(33) 



For M200 > 2 X 10 M0, we find for the F200 versus M200 rela- 
tion a = 1 .682 ± 0.004 and /? = -5 .648 ± 0.002, with a reduced 
of 1.2 for ^ 13,000 clusters divided into 19 bins. This 
slope is slightly steeper than that expected for the self-similar 
case. It is also consistent with the Y-M relations de rived from 
hydrodynamical cluster simulations (e.g . White et al.l l2002t 
Ida Sil va et alTIOOllMotl et alJ2005tlNagai,2006i:lShaw et alJ 
2008). The Y-M relation calculated within a disk of radius 
/?2500 is also consistent with the observed F2500-A^gas(< ^250q1 
relatio n of Bonamente et al., (,2008i) . as shown in .Bode et all 
(l2009h . 
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TABLE 3 

Best- Fit Parameters for the Fa -Ma Relation, with Ma > 2 x 10'"*Mq , 
Fitted Using the Power Law Given in Eouation[33] 
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Fig. 12. — F200 versus M200 relation for all the halos from one octant of this 
simulation with M200 > 2 X IO'^^Mq. The points with error bars represent 
bins of width Alog(M2oo) = 05 and errors on the bin mean. The best-fit 
parameters for this scahng relation and the scatter are given in Tabled 



We also find the scatter in this Y-M relation about the best- 
fit relation. Specifically we calculate 



O'YM ■ 



EZidnYA-lnYA)? 
N-2 



1/2 



(34) 



as discussed in e.g. ^ Shaw et aP (HOO^, where ?a is the pro- 
jected SZ flux measured within a disk of radius R^, 1^ is the 
fitted SZ flux of for a cluster of mass Ma using the best-fit pa- 
rameters for Eq. [33] and N is the total number of clusters. We 
find a scatter for the Y200-M200 relation of about 14%, and a 
slightly larger scatter of about 17% for the Fsoo-^soo relation. 
We find a lower scatter of roughly 13 % for the Y 500-M200 
relation, which has been suggested by Shaw et alJ (12008) to 
have less scatter than the Fioo-A^aoo relation. We list the best- 
fit parameters for Eq. [33] and the scatter for the Y2cn)-M2oo, 
YsoQ-Msoo, and F5()()-M2oo relations in Table |3] 

Concerning the scatter given above, we assume hydrostatic 
equilibrium for the gas; however, we do not assume the same 
entropy level for all the halos. The kinetic energy of the gas 
is derived from that of the dark matter, so for a non-virialized, 
merging halo with a high kinetic energy relative to the cluster 
potential, the gas will also be given a high energy and thus 
temperature. Thus merging of halos is taken into account at 
some level, although a full treatment of the gas physics would 



give more scatter. Also, the repetition of halos as the periodic 
box is tiled to fill the octant, even though these halos are of- 
ten at different redshifts and evolutionary states, may serve to 
slightly decrease the scatter. 

Note that the Y values discussed above are intrinsic to the 
cluster and do not include any contamination from foreground 
or background SZ signals along the line of sight. Inclusion of 
this projection contamination will increase the Y-M scatter 
somewhat, but for clusters with masses greater than roughly 
2 X 10 '4^0, an optimistic mass limit of current SZ surveys, 
this projection contamination should be below ^ 20% over a 
wi de range of possibl e an values (0.7-1.0) and redshifts (z < 
2) (iHolder et al.ll2007t iHallman et al.ll200l . 

3.4. Point Source Contamination of SZ Clusters 

An important question for both determinations of an SZ 
cluster sample's selection function and measurements of the 
SZ power spectrum is the correlation and contamination of 
the SZ cluster signal by infrared and radio galaxies. Observa- 
tions suggest we should expect some preference f or radio and 
infrared galaxies to reside in galaxy cluster s (e.g. ICoble et aLl 
l2007HLin & Mohd2007l;lDaddi et al.12009 '). and these sources 
could potentially fill in the SZ decrement at frequencies below 
the null. This point source contamination of the SZ flux could 
both affect cluster detection and bias a measured SZ flux - 
mass scaling relation. If these effects are significant and unac- 
counted for, then inaccurate cosmological conclusions would 
result from measured SZ cluster samples. Regarding mea- 
surements of the SZ power spectrum and derived values of 
(78, typically a significant number of high flux point sources 
are masked out of microwave maps before determining the 
power spectrum. If these sources reside preferentially in mas- 
sive clusters, a significant portion of the SZ signal could be 
masked out and a biased determination of (t% would arise. 

Here we investigate, via these simulations, the fraction of 
SZ clusters that may be substantially contaminated by radio 
or infrared galaxy populations. We first add the flux of all 
galaxies of a given population residing in a given halo. Then, 
using each galaxy cluster's integrated SZ flux within /?2()0 (in 
units of mJy), we calculate the fraction of halos whose com- 
bined flux from resident galaxies equals 20% or more of the 
halo's SZ flux decrement at 30, 90, and 148 GHz. Such halos 
will be harder to detect by SZ surveys, and their signal is more 
likely to be masked out in attempts to measure the SZ power 
spectrum. We choose a 20% contamination fraction as that is 
roughly where the uncertainty in the expected SZ flux for a 
cluster of a given mass becomes dominated by point source 
contamination as opposed to the intrinsic scatter in the Y-M 
relation or projection effects. 

For the contamination due to radio sources, we find that for 
halos with M200 > 1 x \Q^*Mq, roughly 3% or less of clusters 
with z < 1.6 have their SZ flux decrement at 148 GHz filled 
in by more than 20% (see Fig. [T3T l. At 90 GHz, this fraction 
increases slightly, and less than 4% of clusters with z < 1 .6 
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Fig. 13. — Fraction of halos in this simulation that host radio galaxies 
whose combined flux equals 20% or more of the halo's SZ flux decrement 
at 30 GHz (top), 90 GHz (middle), and 148 GHz (bottom). Bin widths are 
Alog(M2oo) = 0.2, and the error bars assume a Poisson distribution. Red, 
blue, and green dots indicate the redshift ranges of 0.4 — 0.8, 0.8- 1.2, and 
1.2-1.6 respectively. 



have a contamination of 20% or more. At 30 GHz, one may 
expect 20 to 30% of clusters with M2m ^ 1 x lO'^^M© to be so 
contaminated, with a decrease in this fraction as the halo mass 
is increased. These contamination fractions exhibit a weak 
redshift dependence, with slightly larger fractions at higher 
redshifts. 

Infrared galaxies seem to be a significant source of contam- 
ination for SZ clusters at 90 and 148 GHz (see Figs. [T4l and 
[T5] ). However, these simulations suggest that, at 90 GHz, the 
fraction of halos with M2()o > 2.5 x lO'^^M© and z < 1 .2, that 
are contaminated at a level of 20% or more, is less than 20%. 
In all mass bins, the general trend is for the contamination 
fraction to increase at higher redshifts where more infrared 
galaxies reside. There is no appreciable contamination from 
infrared galaxies at 30 GHz. At 148 GHz, less than 20% of 
halos with M201) > 2.5 x 10"*Mq and z < 0.8 have their SZ 
decrements filled in at a level of 50% or more (see Fig.fTSll. It 



should be noted that the above estimates only include galaxies 
that reside in halos and does not include contamination from 
galaxies along a halo's line of sight. In this regard, these frac- 
tions should be viewed as lower limits. The one exception to 
this is for the very high flux infrared galaxies that most likely 
represent a lensed population of galaxies. These galaxies are 
placed in massive halos and are the main source of contami- 
nation in the most massive bins. These sources should more 
properly be placed behind massive halos, however, the con- 
tamination effect in projection will still be similar 

Such contamination fractions as discussed above should to 
be factored into sample completeness estimates and the scat- 
ter in the Y -M relation in order to use SZ detected clusters 
for precision cosmology. The models presented here should 
give a reasonable approximation to the expected levels of con- 
tamination from the different galaxy populations. However, 
these estimates should be cross-checked and made more pre- 
cise with multi- wavelength observations of SZ clusters. 
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Fig. 14. — Fraction of halos in this simulation that host infrared galaxies 
whose combined flux equals 20% or more of the halo's SZ flux decrement at 
90 GHz (top) and 148 GHz (bottom). Bin widths and redshift ranges are the 
same as in Fig. |13| An absence of an indicated fraction indicates no halos in 
that mass and redshift interval. 
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whose combined flux is greater than 10% (green), 50% (blue), or 100% (red) 
of t he h alo's SZ flux decrement at 148 GHz. Bin widths are the same as in 
Fig. 1131 This is shown for the redshift ranges of 1.2— 1.6 (top), 0.8- 1.2 
(middle), and 0.4-0.8 (bottom). 



4. CONCLUSION 

The simulations discussed here have been created for the 
use of the ACT team to test their data analysis pipeline. To 
this end, we have modeled the primary astrophysical signals 
as realistically as possible by matching closely to recent ob- 
servations. We have also expanded the frequency coverage 
beyond the ACT bands and created maps and catalogs in an 
easily accessible format. As such, other current and upcom- 
ing microwave background experiments may find these simu- 
lations useful for a variety of purposes. 

From the gas prescription in these simulations, which has 
been matched to recent X-ray observations, we find a scaling 
relation between Compton y parameter and mass that is only 
slightly steeper than self-similar and has a scatter of 14% to 
17% around the best-fit relation for Fioo-A^aoo and Isoo-Msoo 
respectively. This low intrinsic scatter suggests that measure- 
ment of the cluster SZ flux may be a relatively clean clus- 
ter mass proxy. The correlation of radio galaxies with SZ 



clusters suggests that at 148 GHz (90 GHz), for clusters with 
M20() > IO'^^Mq, less than 3% (4%) of these clusters will have 
their SZ decrement filled in by 20% or more. We find the con- 
tamination levels higher for infrared galaxies. However, at 90 
GHz, less than 20% of clusters with M200 > 2.5 x IO'^^Mq 
and z < 1 .2 have their SZ decrements filled in at a level of 
20% or more. At 148 GHz, less than 20% of clusters with 
M20() > 2.5 X IO^'^Mq and z < 0.8 have their SZ decrements 
filled in at a level of 50% or larger. We also find that a pop- 
ulation of very high flux infrared galaxies, which may have 
been lensed, contribute most to the SZ contamination of very 
massive clusters at 90 and 148 GHz. Multi-wavelength obser- 
vations of SZ clusters will help to make these contamination 
estimates more precise. This information is crucial for deter- 
mining SZ cluster selection functions and the actual scatter in 
the observed Y-M relation, which in turn are vital to efforts 
to extract cosmology from SZ cluster counts. The inclusion 
of the lensing of the microwave background in a manner in- 
ternally consistent with the large and small-scale structure is 
necessary to verify the accurate recovery of cosmology from 
cross-correlation studies between the lensed microwave back- 
ground and tracers of large-scale structure. Moreover, the 
simulations as a whole, matched to the most recent observa- 
tions, will inform, through a variety of means, measurement 
of the angular power spectrum, higher point correlation func- 
tions, and the cosmology thus derived. 

Cosmic microwave background experiments are entering an 
exciting time of finer resolution and higher sensitivity which 
will open up the Sunyaev-Zel'dovich Universe, microwave 
background lensing, an assortment of cross-correlation stud- 
ies, and deeper probes of inflation through tighter parameter 
constraints. These simulations should provide a useful tool to 
aid analyses in these new and promising areas of discovery. 
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